Exchange and correlation as a functional of the local density of states 
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A functional Exc[p{r, e)] is presented, in which the exchange and correlation energy of an electron 
gas depends on the local density of occupied states. A simple local parametrization scheme is 
proposed, entirely from first principles, based on the decomposition of the exchange-correlation 
hole in scattering states of different relative energies. In its practical Kohn-Sham-like form, the 
single-electron orbitals become the independent variables, and an explicit formula for the functional 
derivative is obtained. 
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I. INTRODUCCTION 

Density functional theory^'^ (DFT) has become an es- 
sential tool in the study and simulation of molecules 
and condensed systems. Its evolution from the initial 
local density approximation (LDA) to the generalized 
gradient approximation"^"^ (GGA), the mixing of exact 
exchange^'^ and, more recently, to the meta-GGAs^'^ has 
brought a substantial increase in accuracy. However, this 
evolution is not systematic and the holy grail of "chemi- 
cal accuracy" (~ 1 Kcal/mole for reaction energies) still 
appears rather distant. Perhaps more importantly, basic 
difficulties remain unsolved, specially for strongly corre- 
lated and magnetic systems, as well as for excited states. 
Among the latter, the systematic underestimation of the 
insulating band gap is particularly troublesome for the 
simulation of solids. 

In this work, I propose a functional for the exchange 
and correlation energy E,j:c, which might potentially 
overcome some of these difficulties. Like in the meta- 
GGAs, its basic ingredient is not just the electron den- 
sity. Rather, E^c depends on a richer and more informa- 
tive function, namely the local density of occupied states. 
Since this function changes much more drastically and 
transparently than the total density for excited states, 
the functional can in principle describe these in addi- 
tion to the ground state. Furthermore, the functional 
form can be rationalized naturally in terms of electron- 
electron scattering states, what in principle allows also 
to parametrize the functional in an entirely ab initio and 
unambiguous way. 



II. DENSITY OF STATES FUNCTIONAL 

Let us first define the local density of occupied kinetic- 
energy states as 
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p{r, e) = ^ J d'^r rjv)P S'{r, - r) 



with 
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where means real part and I will omit spin degrees of 
freedom to simplify the notation. Notice that the first 
two moments of p(r, e), 
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are the electron and the kinetic energy densities. The 
former is the functional variable in all DFT's and the 
latter is used as an additional auxiliary function in meta- 
GGAs*'^ (with p(r) and r(r) obtained from a noninter- 
acting wave function, see below). These may thus be con- 
sidered as particular cases of a more general functional 
E[p{r,e)]. 

Notice that, since the wavefunction ^'(ri, tat) is de- 
termined^ by the electron density p(r), the density of 
states p{r,e) is also fully determined by p(r), and there- 
fore it formally does not add any variational freedom. 
However, in practice p(r, e) provides explicitly much more 
information on the electronic state than p(r), and there- 
fore it allows for a more accurate parametrization of 
the exchange and correlation energies. On this respect, 
p{r, e) is a reasonable half way in complexity between 
p{r) and \l/(ri, rjv). Thus, I express the total energy 
as 

E[p{r,e)] = T[p{r,e)] + UH[p{r)]+E.,4p{r,e)] 

dh p{r) Vextir) (5) 
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where T is the kinetic energy 

r[p(r,e)] = I d^rrir), 
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Uh is the Hartree energy 

aMMr)l = //cl«r.Vf!« (7) 

Vexti^) is the external potential, and Exc[p{y^,e)] is some 
approximation to the exchange and correlation energy, 
which includes everything else. 

Although Eq. (6) is in principle the total kinetic en- 
ergy, in practice we will use a density p(r, e) obtained 
from a single-determinant wave function and there- 
fore Exc[p{j^, e)] must also contain the difference in kinetic 
energy between and the exact wave function. More 
precisely, let us introduce another functional 

E°[p°{r,e)] = T[pO(r,e)] + [/H[p(r)]+i;L[P°(r,e)] 

+ j d^r p(r) Vext(r) (8) 

where 

p%r,e) = f2\Mr)\'S{el{r)-e), (9) 

a=l 



The advantage of p(r, e) as a functional variable is that 
it provides much more explicit information on the system 
state than p(r). Thus, it naturally leads to an energy and 
orbital-dependent exchange-correlation potential, much 
like quasiparticle self energies, thus offering hope of im- 
proving the description of excited states and band gaps. 
Also, since the density of states is very different in local- 
ized and extended systems, it offers the possibility to ad- 
dress electron localization and self-interaction problems. 
In a different context, the relevance of spectral functions 
in the description of electron localization and strong cor- 
relation has been stressed by dynamical mean field thc- 
oryio. Finally, the success of meta-GGA functionals in 
improving atomization energies gives us some confidence 
that the proposed functional may also improve consider- 
ably on formation and reaction energies. On the other 
hand, it might be suspected that, since the local density 
of states is much richer in structure than the total den- 
sity, a functional dependent on it should be also much 
harder to parametrize. I will show that, at least in prin- 
ciple, this is not the case, and that unambiguous and 
physically rooted parametrizations of E^c are possible. 
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V'a(r) 
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and V'a(r) are single-electron, orthonormal, occupied or- 
bitals. 

For a given density p{r), the functions /o(r, e) and 
p'^(r,e) are different. Furthermore, since the same den- 
sity p{r) can result from different sets of V'a(r)'s, and 
therefore from different functions p^{r, e), there is an in- 
trinsic indetermination in the functional E'^ . An example 
of such indetermination is a rotation in the space of occu- 
pied orbitals, which leaves invariant p(r) but not /C°(r, e). 
The indetermination can be removed by specifying an ex- 
plicit way of obtaining ipa{i^) from p(r). A natural choice 
is to minimize T[p'^{r, e)], under the constraint (3). This 
leads to the Kohn-Sham equations 



(11) 
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-V^ + v(r)-eaU'a(r)=0 



where v{r) is the Lagrange multiplier associated to the 
constraint (3), and £„ is that associated to the normal- 
ization constraint for ipa- However, I will simply consider 
the indetermination in /3°(r, e) as an additional freedom 
to achieve a practical and accurate parametrization of 
E°[p°{r,c)]. 

In Eqs. (5) and (8), E[p{r,€)] and E°[p°{r,€)] are for- 
mally assumed to be infinite for all functions, p{r, e) and 
/9°(r, e), which cannot be expressed in the forms (1) and 
(9), respectively. In practice, this condition needs not 
be dealt with, since we will use the standard Kohn-Sham 
procedure^ of minimizing the total energy as a function of 
V'a(r) rather than of /9°(r, e). In what follows, I will refer 
only to E^^, and p"{r, e), but I will drop the superscript 
to simplify the notation. 



III. FUNCTIONAL FORM 
Here I propose the local functional form^^ 

d^r II de de' K^r), e, e') e) e') 

(12) 

where the kernel K^c is an average of the interaction 

energy between two electrons with wave vectors k of fixed 
size k = y/2mee/h but random direction k: 

= IJ dqexc{p,q) Piq\e,e') (13) 

where the factor 1/2 corrects double counting and 
P{q\e,e') is the probability that the two electrons, with 
kinetic energies e = h^k'^/2me and e' = h^k'"^ /2me, have 
relative momentum hq = h\k — k' \ / 2: 



, , d'k d^k' f \k-k'\ 

2q/{kk') if \k - k'\ <2q<\k + k'L^. 
otherwise. ^ ' 



&xc{p, q) is the energy change of an electron pair, with 
relative momentum hq, in a homogeneous gas of den- 
sity p (assumed unpolarized for simplicity), when their 
antisymmetry is imposed and their interaction potential 
e^/r is switched on. It may be written as e^dp, q) = 
^x{q) + ec(p, q)- Expanding the plane waves in spherical 
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harmonics 



°° roo 2 

eM = V(2I + 1) {si - 1) / dr 4TTr^- jf{qr) 

Jo I" 
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■ V 



(15) 



where the factor s; = 1 — (— 1)'/2 accounts for the an- 
tisymmetry of equal-spin electron-pair wavefunctions^^. 
Cc may be approximated by 

OO POO 

ec{p, q) = ^(2/ + 1) SI / dr Anr^ 

1=0 -^0 
X [(t>i{p,Q,r) Hi (t)i{p,q,r) - ji{qr) Hi ji{qr)], {16) 

where Hi =Ti + e^/r and 



Sr^ 2/ir2 ' 



(17) 



with ^ — me/2. <j)i{p, q, r) and ji{qr) are the interacting 
and noninteracting scattering wavefunctions: 



Ti + Vint{p, r) 



2ii 



Mp,q,r) = o, (18) 



and ji{qr) is the solution of the same equation with 
Vint = 0, i.e. a spherical Bessel function. Vint{p,'r) is an 
effective interaction potential for the screened Coulomb 

repulsion, designed to reproduce accurately the pair cor- 
relation function of the homogeneous electron gas^"^: 

g{p,r)=y2{2l + l) SI dqP{q\p)<i>f{p,q,r), (19) 

where P{q\p) is the probability of finding an electron pair 
with relative momentum ?iq in a homogeneous gas of den- 
sity p. For qKkp = (37r^p)^/^: 

P{q\p) = j de de' P{q\e,e') P{e\p) P{e'\p) 



0^ 

tvp tVjp tvp 



(20) 



where P{e\p) = 3k/kp is proportional to the free-electron 
density of states. Eqs. (18-19) define a problem which 
is presently under very active rcscarch^^'^"'"^'^. Fur- 
thermore, there a other possible alternatives for Cc in 
Eq. (16), like using e^/r instead of Hi and integrating 
in the coupling constant of the electron-electron interac- 
tion. Therefore, the actual parametrization of Vint{p,r) 
and exc{p, (?) will be left for a future work. 

A problem with the above approach is that, while the 
local kinetic energies ea(r) are negative in classically- 
forbidden regions, they are always positive in a homo- 
geneous gas. This problem may be overcome by using 
the Airy gas proposed by Kohn and Mattsson-^° in which 



the total one-electron effective potential is a constant- 
force ramp v{r) = i ■ (r — Tq). Since a constant force 
does not affect the relative motion of two equal particles, 
it is a good approximation to use the same Cxc of the 
homogeneous gas and to modify only the probability P 
of equation (13) for the Airy gas. To do this, we can 
use the one-to-one relationship^" (p(r), V/!)(r)) <^ (f, ro), 
making P and K^c dependent on |Vp| in addition to p. 
In fact, this procedure can be used also to derive GGAs 
and meta-GGAs from exc{p,q)- 



IV. FUNCTIONAL DERIVATIVE 

To minimize the total energy or to perform ab initio 
molecular dynamics, we need the functional derivative 
SExc/Sip^{r)'^^''^'^ . To calculate it, it is convenient to as- 
sume that the occupied orbitals tpa{r), density p(r), and 
density of states p{r, e) are given in a mesh of points 
and discrete energies e. Their gradients and laplacians 
can then be calculated as 



Vpi = Vp(ri) = 



ijPj 



(21) 



(22) 



where pi = p{ri),tpja = tpa{i'^j), and gij,Lij are coeffi- 
cients of some unspecified finite-difference formulas. We 
can then replace the functional derivatives by ordinary 
partial derivatives and recover the continuous limit at 
the end. Using atomic units {Ti = e = rrie = 1) and omit- 
ting for brevity the grid increments, Ax'^ and Ae, when 
replacing the integrals J d^r and / de by grid sums, the 
discretized equations are 



3 ee' 



(23) 



where K^^^' = Kxc{p{vj), |Vp(rj)|, e, e'), 

Pje = ^^ja^3a^jae, (24) 
a 

with Sjae a discretized delta function (remember that e 
are now discrete energies) like 

^ if \eja — e| < Ae 
otherwise, 

(25) 

which must verify that J2e ~ ^ ^^ja , where 



6jae = S{eja-e) 



Ae I ^ Ad 



2i>jo 



I will now consider formally ij)*-^ and V'ja as independent 
variables to perform the derivatives 



dP] 



* — ^ij'4'jcft 



(27) 
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where Sij is the usual Kronecker delta, 

-^ia, (28) 



= dij^pjadjae + 



+ 



4 de., 



and 



+2 5I^^«'^J^ 



(31) 



Substituting (27), (28), and (30) into (31) we obtain 
dE^c I 



2 ^ SiaeVi^ Ipi, 



dSi, 



2V'* 



Wear 



l^Uj^j^Y^^vj.. (32) 

1 p 3^ 



I have defined 



Vie = ^ Kiee'Pie' 



(33) 
(34) 



where I have used that Lj, = Lij and g^i = —gij. Re- 
turning to the continuous limit, 

= (^«(r) + 2 y de ^(e„(r) - e) ^r, e)) ^a(r) 

^„(r)W*(r) /■ d,?(e^(r)-e) 
+ 2V>S(r) 7 - "^"'^^ 



de„(r) 

_d5(ea(r) - e) 



w(r) 



^(r,e) (35) 



dede' 



dp{v) 



P(r, e)p(r, e') 



-V 



/ Vp(r) 
V|Vp(r)| 

' ,aii'xc(r,e,e') 



p(r,e)p(r,e') , (36) 



w(r, e) 



Kxc{r, e, e') p(r, e' 



(37) 



where I am using the shorthand notation K^d'^. e, e') = 
-?'^a;c(p(i'), |Vp(r)|, e, e'). Using integration by parts in 
(35) to integrate out the delta functions 



(«(r) + 2u(r,e„(r)))Va(r) 

Vc«(r)V2V';(r) dv{r,eair)) 



+ 



2V*(r) 



V'a(r) 



5e„(r) 

dv{r,€a{r)) 
dea{r) 



(38) 



In practice, v{r,ea{r)) and its derivative will be interpo- 
lated from the values at a number of discrete energies, i.e. 
v{ri,ea{i^i)) = WieCe,(ea(rj)), where the interpolation 
functions arc determined by, e.g. Lagrange or spline 
interpolation formulas. Comparison with Eq. (32) shows 
that the discrete delta values Siae, used to distribute the 
contribution of state a to p{ri,e), at the discrete ener- 
gies e, are simply 5iae = Ce(ea(rj)). Thus, for example, 
Eq. (25) corresponds to a simple linear interpolation for 
v(ri,ea{ri)). Though Eqs. (36-38) arc considerably more 
complicated than the functional derivatives in LDA or 
GGA functionals, the gradients and laplacians can be 
calculated efficiently using fast Fourier transforms or fi- 
nite differences, at a quite affordable computational cost. 



V. CONCLUSIONS 

In conclusion, I have proposed a functional form for 
the exchange and correlation energy, with the potential 
to describe more accurately excitation energies and band 
gaps than the present forms. The basic ideas, embodied 
in Eq. (12), are that the interaction between two electrons 
depends on their relative motion, and that the exchange 
and correlation energy of an electron depends on its ki- 
netic energy. These ideas are physically transparent and 
intuitively appealing, and therefore they may be expected 
to allow for significant improvements in the description 
of the inhomogeneous electron gas. 
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